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We examine the network of forces to be expected in a static 
assembly of hard, frictionless spherical beads of random sizes, 
such as a colloidal glass. Such an assembly is minimally con- 
nected: the ratio of constraint equations to contact forces 
approaches unity for a large assembly. However, the bead 
positions in a finite subregion of the assembly are underde- 
termined. Thus to maintain equilibrium, half of the exterior 
contact forces are determined by the other half. We argue that 
the transmission of force may be regarded as unidirectional, 
in contrast to the transmission of force in an elastic mate- 
rial. Specializing to sequentially deposited beads, we show 
that forces on a given buried bead can be uniquely speci- 
fied in terms of forces involving more recently added beads. 
We derive equations for the transmission of stress averaged 
over scales much larger than a single bead. This derivation 
requires the Ansatz that statistical fluctuations of the forces 
are independent of fluctuations of the contact geometry. Un- 
der this Ansatz, the d(d + 1) /2-component stress field can be 
expressed in terms of a d-component vector field. The proce- 
dure may be generalized to non-sequential packings. In two 
dimensions, the stress propagates according to a wave equa- 
tion, as postulated in recent work elsewhere. We demonstrate 
similar wave-like propagation in higher dimensions, assuming 
that the packing geometry has uniaxial symmetry. In macro- 
scopic granular materials we argue that our approach may be 
useful even though grains have friction and are not packed 
sequentially. 

PACS numbers: 46.10.+Z, 83.70.Fn 



I. INTRODUCTION 

The nature of the forces within a static pile of grains 
has proven more subtle than one might expect. Such a 
pile is an assembly of many hard, spheroidal bodies that 
maintain their positions via a balance of gravitational 
forces and contact forces with their neig hbors @,§. On 
the one hand, determining these forces is a prosaic equi- 
librium problem. Since the number of grains is large, 
the long-established notions of continuum solid mechan- 
ics appear applicable. On the other hand, a pile of grains 
or beads is not a solid. The forces between beads are 
more problematic than those between the atoms of a con- 
ventional solid. These latter forces are smoothly varying 
on the scale of the separation and they arise from a po- 
tential energy that includes attraction. The forces on a 
grain are different. First, they vary sharply with inter- 
particle distance, and there is no attraction. Second, the 
frictional part of a contact force is not determined by 
the macroscopic positions of the grains. Rather, it de- 



pends on the how each contact was formed. The resulting 
macroscopic behavior of the pile is also clearly different 
from that of a conventional solid. Arbitrarily slight forces 
can disturb the pile, so that the notion of stable equilib- 
rium is suspect. Despite these complexities, we expect 
the mechanics of a granular pile to be universal. Hard, 
round grains appear to form piles of the same nature in- 
dependent of their composition or detailed shape. We are 
led to think of these as nondeformable objects that exert 
normal forces of constraint and transverse forces limited 
by Coulomb's static friction limit. 

Recently, a puzzling discovery has underlined the sub- 
tlety of the forces in a conical heap of poured sand || . 
The supporting force under the center, where the pile is 
deepest, is not maximal. Instead, the maximal force oc- 
curs along a circle lying between the edge and the center 
of the pile. From this circle the force decreases to a min- 
imum at the center. In order to explain this puzzling 
"central dip" , a number of inventive approaches have 
been taken. Some PJ seek to account for the minimum 
qualitatively by viewing the pile as a stack of concentric 
"wigwams" , whose sloping sides support the load. Others 
p] have shown that the central dip is compatible with the 
conventional continuum mechanics, in which the pile is 
viewed as a central elastic zone flanked by an outer plas- 
tic zone which is at the Coulomb static friction limit. A 
third group (|,[l4|,[l5|] has argued that granular material 
requires a new constitutive law, a homogeneous, local, 
linear constraint on the stress arising from the packing 
geometry. We shall call it a "null stress" law. Their pro- 
posed law gives a continuum mechanics as simple as that 
of a liquid or a solid, yet different from either. The hall- 
mark of this law is the hyperbolic equations governing the 
transmission of forces. Hyperbolic equations, such as the 
wave equation, obey causality. The wave at the present 
is unaffected by the future. In the null-stress picture, 
the vertical direction plays the role of time. Accordingly, 
forces at a given point in the pile are only influenced by 
forces above that point. Force propagates as in a travel- 
ing wave. The transmission is unidirectional, in contrast 
with conventional continuum elasticity. The equations of 
continuum mechanics are elliptical. According to these, 
a contact within the pile should be influenced by all the 
forces above or below it, as sketched in Figure El 

A separate approach has given indirect evidence for the 
unidirectional transmission of forces. Coppersmith's 
heuristic q model aims to account for the point-to-point 
variability of the contact forces. It imposes a unidirec- 
tional prescription for determining the downward forces 
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from one grain in terms of the forces acting on it from 



above. Both this model and refinements |tL8|J19 of it yield 
an exponential falloff of probability for large forces. This 
exponential falloff agrees well with measured distribu- 
tions This exponential falloff contrasts with 
the Gaussian falloff expected for a heterogeneous elastic 
solid. 

In this paper we grapple with the relationship between 
the conventional elastic view and the newer null-stress 
picture. Our work builds on several recent studies of the 
relationship between the connectivity of a structure and 
the force transmission in it. Alexander's recent review 
p2j has explored the nature of unconstrained degrees of 
freedom in a minimally connected or isostatic pi] net- 
work. Ball and Edwards jll) have explored force trans- 
mission through minimally-connected networks, assum- 
ing a fixed co-ordination number for all particles. They 
have shown that such lattices can have constitutive equa- 
tions of the null-stress form. Our main aim here is to 
broaden the class of systems that must show unidirec- 
tional force transmission, as required by the null-stress 
picture. 

We present the discussion of the problem of stress 
transmission in a granular packing on several levels of 
generality. In sections P |fv| we focus on properties of the 
system of frictionless, spherical beads. Possible experi- 
mental realizations are hard-sphere colloidal dispersions 
f24| , ^5f or weakly-deformed droplet emulsions |26| . First, 
in Section [il] we use general counting arguments to show 
that such a packing is minimally coupled. We then re- 
late this fact to the inadequacy of elastic description for 
such a system. For a small subsystem within the pile a 
counting of equations and unknowns shows that approx- 
imately half of the surface forces transmitted from out- 
side the subsystem are redundant. In equilibrium these 
cannot be independently specified, but have a fixed re- 
lation to the other half of the forces. This requirement 
for balance amongst many forces implies the existence of 
soft modes — infinitesimal deformations with no restoring 
force. In the continuum limit, the soft modes impose 
conditions on the stress field of the null-stress form. 

In order to obtain the particular form of such macro- 
scopic description, we limit further discussion to so-called 
sequential packing specified in Section III. In Section 



IV we give a microscopic prescription for determining 
the contact forces in a unidirectional way, and develop a 
green's-function formalism for determining the forces in 
the d-dimensional sequential packing. This picture al- 
lows one to decouple the geometric features of the pack- 
ing from the pattern of transmitted forces. In Section IV 
we also explore the macroscopic consequences of these 
force laws, leading to an expression for the stress tensor 
with d variables rather than the d(d + I)/2 variables of 
a general stress tensor. In two dimensions, our formal- 
ism places a constraint on the stress, whose form agrees 
with the null stress law of Wittmer et al. || . In higher 
dimensions, the constraints on the stress also lead to a 
unidirectional equation for the transmission of stress in 



the form of a wave equation. In Section 6 we consider the 
relevance of our findings for real granular piles, which are 
not sequentially packed and which have friction. Friction 
can alter the transmission of forces qualitatively, restor- 
ing the elastic transmission of stress. 



II. UNDERDETERMINATION WITHIN A 
FRICTIONLESS PACK 

In this section we consider how the constraints inherent 
in the packing of impenetrable, frictionless beads deter- 
mine the contact forces between the beads. Forces in fric- 
tionless packs have been studied by simulation . 
Theoretical properties of the forces have been established 
for simplified systems jll]-[l3| . We begin with a summary 
of the well-recognized enumeration of equations and un- 
knowns, as discussed, eg. by Alexander [^2| and by the 
Cambridge group |Tl]] . We then consider the role of ex- 
ternal forces acting on a subregion of the pack, and show 
that a subset of these suffices to determine the others. 

For definiteness we consider a system of rigid spherical 
beads whose sizes are chosen randomly from some contin- 
uous distribution. We suppose that the diameters of all 
the spheres are specified and the topology of their pack- 
ing is given. That is, all the pairs of contacting beads 
are specified. The topology can be characterized by the 
average coordination number, i.e.,the average number of 
nearest neighbors, Z = 2N C /M (here N c is the total num- 
ber of contacts, and M is the number of beads). The 
necessary condition for the packing of a given topology 
to be realizable in a space of dimensionality d is that the 
coordinates of the bead centers, x Q satisfy the following 
equations, one for each pair of contacting beads a and (3: 



(x Q - = (R a + R/s 



(1) 



Here R a , Rp are the radii of the beads. There are N c such 
equations (one for each contact) , and Md variables (d for 
each bead). The number of equations should not exceed 
the number of variables; otherwise, the co-ordinates x Q 
are overdetermined. Thus, the coordination number of 
a geometrically-realizable packing should not exceed the 
critical value of 2d: Z < 2d. We assume all the equations 
imposed by the topological constraints to be indepen- 
dent. If they were not independent, they would become 
so upon infinitesimal variation of the bead sizes. For in- 
stance, the hexagonal packing in two dimensions has the 
coordination number 6 which is higher then the critical 
value, (Z) — 4; but the extra contacts are eliminated 
by an infinitesimal variation of the bead diameters. In 
other words, the creation of a contact network with coor- 
dination number higher than 2d occurs with probability 
zero in an ensemble of spheres with a continuous distri- 
bution of diameters. We shall ignore such zero-measure 
situations henceforth. 

The above consideration gives the upper limit on the 
average coordination number, Z. The lower limit can be 
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obtained from the analysis of mechanical stability of the 
packing: it gives a complementary inequality: Z > 2d. 
We will consider a packing to be mechanically stable if 
there is a non-zero measure set of external forces which 
can be balanced by inter-bead ones. The packing of fric- 
tionless spheres is always characterized by (Z) — 2d , as 
we now show. Stability requires that the net force on each 
bead be zero; there are Md such equations. The forces 
in these Md equations are the iV c contact forces. The 
Md equilibrium conditions determine the magnitudes of 
the N c contact forces. (Their directions are determined 
by the geometry of the packing.) The number of equilib- 
rium equations Md should not exceed the number of force 
variables N c ; otherwise these forces would be overdeter- 
mined. Thus Md < N c , or ~Z > 2d. To avoid both 
overdetermined co-ordinates and overdetermined forces, 
we must thus have Z = 2d. 

Similar counting arguments have been discussed pre- 
viously H^JlJI . A subset of them have been applied to 
granular packs with friction [jyj. Here we emphasize a 
further feature of a frictionless bead pack that has not 
been well appreciated: the co-ordinates and forces within 
a subregion of a large bead pack are necessarily under- 
determined. Quantifying this indeterminacy will play an 
important role in our reasoning below. To exhibit the 
indeterminacy, we consider some compact region within 
the packing, containing M' beads. This unbiased selec- 
tion of beads must have the same average co-ordination 
number Z as the system as a whole: Z = 2d. Let iV ex t 
be the number of contacts of this sub-system with ex- 
ternal beads, and Ni nt be the number of the internal 
contacts. The average coordination number Z can be 
expressed Z — (N cxt + 2Ni nt )/M' (any internal contact 
should be counted twice). Since there are M'd equa- 
tions of force balance for these beads, one is able to 
determine all N cxt + N m t contact forces in the system, 
whenever M'd = N cxt + N[ nt . Evidently, if the forces on 
the iVext contacts are not specified, the internal forces 
cannot be computed: the system is underdetermined. 
The number of external forces Nq required is given by 
No = M'd — Nix& . This No may be related to the average 
co-ordination number Z : 



N n = M' 



d- Z - 
2 



N c , 



(2) 



We now observe that the quantity in [...] vanishes on 

average. This is because the average of Z for any subset 
of particles is the same as the overall average. There is no 
systematic change of Z with M' . Thus if one half (on av- 
erage) of mutually-independent external forces is known 
(let us call them "incoming" ones), the analysis of force 
balance in the region enables one to determine all the re- 
maining forces, including the other half of external ones 
("outgoing"). We are free to choose the incoming con- 
tacts at will, provided these give independent constraint 
equations. 



This observation supports the unidirectional, propa- 
gating stress picture, discussed in the Introduction. In- 
deed, one can apply the above arguments to the slabs 
of the packing created by cutting it with horizontal sur- 
faces. In a given slab of material, we choose the forces 
from the slab above as the incoming forces. According to 
the preceding argument, these should determine the out- 
going forces transmitted to the slab beneath. This must 
be true provided that the constraints from the upper slab 
are independent. 

Such force transmission contrasts with that of a solid 
body, as emphasized in the Introduction. If a given set 
of forces is applied to the top of a slab of elastic solid, 
the forces on the bottom are free to vary, provided the 
total force and torque on the slab are zero. Yet in our 
bead pack, which appears equally solid, we have just con- 
cluded that stability conditions determine all the bottom 
forces individually. In deducing this peculiar behavior, 
we did not exclude tensile forces; we may replace all the 
contacts by stiff springs that can exert strong positive or 
negative force, without altering our reasoning. In this 
sense our result is different from the recent similar result 
of Moukarzel E^]. The origin of the peculiar behavior 
lies in the minimal connectivity of the beads. 

In a subregion of the minimal network, the constraints 
can be satisfied with no internal forces. Moreover, numer- 
ous (roughly iV 0Xt /2) small external displacements can be 
applied to the subregion without generating proportional 
restoring forces. We call these motions with no restoring 
force "soft modes" . If we replace these external displace- 
ments with external forces and require no motion, com- 
pensating forces must be applied elsewhere to prevent 
motion of the soft modes. If the applied forces perturb 
all the soft modes, there must be one compensating force 
for each applied force to prevent them from moving — on 
average N cxt /2 of them. The subregion is "transparent" 
to external forces, allowing them to propagate individu- 
ally through the region. 

This transparent behavior would be lost if further 
springs were added to the minimal network, increasing 
Z. Then the forces on a subregion would be determined 
even without external contacts. The addition of exter- 
nal displacements would deform the springs, and pro- 
duce proportional restoring forces. There would be no 
soft modes, and no transparency to external forces. 

A simple square lattice of springs provides a concrete 
example of the multiple soft modes predicted above. Its 
elastic energy has the form 



H = K dx&y [{u xx f + (u vv f 



(3) 



This functional does not depend on u xy , thus there are 



shear deformations (u x 



,mj 



0) which cost no elas- 



tic energy. This means that the stress field should be 
orthogonal to any such mode, i.e., 







(4) 
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where u a xx — u yy = 0, and u xy is an arbitrary func- 
tion of (x; y). The above equation implies that a xy = 0, 
j.e.,the principal axes of the stress tensor are fixed and 
directed along x and y. This provides a necessary closure 
for the standard macroscopic equation of force balance, 

= fL (5) 

here f ext is an external force. Since a xy = the two 
unknown components of the stress field, a xx and a yy 
propagate independently along the corresponding char- 
acteristics, x — const and y = const: 

d x a xx = /* xt (6) 

d«o-yy = f^ t (?) 

The propagation of the solution along characteristics is 
a property of hyperbolic problems such as wave equation. 
The above equations without external force imply that 
each component of the stress tensor a satisfies a wave 
equation of the form 

where t = x + y and s = x — y. Thus, the fact that 
the original elastic energy has the soft modes results in 
hyperbolic, rather than elliptic equations for the stress 
field. One has now to specify the surface forces (or dis- 
placements) at a single non-characteristic surface — a line 
not parallel to x or y — in order to determine the stress 
field in all the sample. 

A frictionless granular packing behaves like this exam- 
ple: they both are minimally coupled; they both have 
soft modes; they both have unidirectional propagation. 
In both examples only the surface of the sample stabi- 
lizes the soft modes. The above consideration of regular 
lattice can be easily extended to the case of arbitrary 
angle between the characteristic directions, x and y. In- 
stead of starting with a square lattice, we could have 
applied a uniform x — y shear, altering the angle between 
the horizontal and vertical springs. The reasoning above 
holds for this lattice just as for the original square one. 

The nature of the soft modes in a disordered bead pack 
is less obvious than in this lattice example. We have not 
proven, for instance, that all the forces acting on the top 
of a slab correspond to independent soft modes, which 
determine the forces at the bottom. Otherwise stated, 
we have not shown that the soft modes seen in the mi- 
croscopic displacements have continuum counterparts in 
the displacement field of the region. However, the follow- 
ing construction, like the lattice example above, suggests 
that the soft modes survive in the continuum limit 

To construct the pack, we place circular disks one at a 
time into a two-dimensional vertical channel of width L. 
(Such sequential packings will figure prominently in the 
next section.) Since the disks are of different sizes, the 



packing will be disordered. We place each successive disk 
at the lowest available point until the packed disks have 
reached a height of order L, as shown in Figure |2|. We 
now construct a second packing, starting from a channel 
of slightly greater width L + 6. We reproduce the packing 
constructed in the first channel as far as possible. We 
use an identical sequence of disks and place each at the 
lowest point, as before. There must be a nonvanishing 
range of S for which the contact topology is identical. 
The motion of the wall over this range is thus a soft 
mode. As the side wall moves, the top surface will move 
by some amount e, proportional to S. Now, instead of 
holding the side wall fixed, we exert a given force f x on it. 
Likewise, we place a lid on the top, remove gravity, and 
exert another force f y . Evidently unless f x /f y = e/5, 
a motion of the soft mode would result in work, and 
the system would move. Thus f y plus the condition of 
no motion determines f x . This condition translates into 
an imposed proportionality between the stresses a yy and 
a xx , as in the lattice example above. The soft modes 
have continuum consequences. 

III. SEQUENTIAL PACKING UNDER GRAVITY 

In the previous section we have shown that a packing of 
frictionless spherical beads is an anomalous solid from the 
point of view of classical elastic theory. The fact that the 
average coordination number in such a packing is exactly 
2d for the infinite system supports unidirectional, prop- 
agating stress. Now we elaborate this concept in more 
detail, by deriving particular laws for microscopic and 
macroscopic force transfer adopting a particular packing 
procedure. We suppose that the beads are deposited one 
by one in the presence of gravity. The weight of any 
new sphere added to the existing packing must be bal- 
anced by the reactions of the supporting beads. This is 
possible only if the number of such supporting contacts 
is equal to the dimensionality d. Any larger number of 
contacts requires a specific relationship between the sizes 
and coordinates of the supporting beads, and and thus 
occurs with vanishing probability. As a result, the even- 
tual packing has an average coordination number 2d, like 
any stable, frictionless pack. In addition, it has a further 
property: a partial time-like ordering. Namely, among 
any two contacting beads there is always one which has 
found its place earlier than the other (the supporting 
one), and any bead has exactly d such supporting neigh- 
bors. Note that the supporting bead is not necessar- 
ily situated below the supported one in the geometrical 
sense. The discussed ordering is topological rather than 
spatial. 

One could expect that although any bead has exactly d 
supporters at the moment of deposition, this may change 
later. Specifically, adding an extra bead to the packing 
may result in the violation of positivity of some contact 
force in the bulk [jlij. This will lead to a rearrange- 
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merit of the network. For the moment we assume that 
the topology of the sequential packing is preserved in the 
final state of the system, and return to the effect of rear- 
rangements in Section [V|. 

The partial ordering of the sequential packing consid- 
erably simplifies the calculation of the force distribution. 
Indeed, any force applied to a bead can be uniquely de- 
composed into the d forces on the supporting contacts. 
This means that the force balance enables us to deter- 
mine all the "outcoming" (downward) forces if the "in- 
coming" ones are known. Therefore, there is a simple uni- 
directional procedure of determination of all the forces in 
the system. Below, we use this observation to construct 
a theory of stress propagation on the macroscopic scale. 

IV. MEAN-FIELD STRESS 

We will characterize any inter-bead contact in a se- 
quential packing with a unit vector directed from the 
center of supported bead a toward the supporting one /?, 

X /3 X Q / n x 

n Q /3 = 1 (9) 

I X/3 -x Q | 

The stress distribution in the frictionless packing is given 
if a non-negative magnitude of the force acting along any 
of the above contact unit vector is specified. We denote 
such scalar contact force as f a p 

The total force to be transmitted from some bead a to 
its supporting neighbors is the sum of all the incoming 
and external (e.g. gravitational) forces: 

F Q = (fcxt)a + ^ npaffia (10) 

/3(->«) 

Here j3(— > a) denotes all the beads supported by a. 
Since there are exactly d supporting contacts for any bead 
in a sequential packing, the above force can be uniquely 
decomposed onto the corresponding d components, di- 
rected along the outcoming vectors n a7 . This gives the 
values of the outcoming forces. The /'s may be com- 
pactly expressed in terms of a generalized scalar product 
<-l-> Q : 

f ai = (F Q |n Q7 ) Q (11) 

The scalar product (...|...) a is defined such that 
(n Q7 |n Q7 /) = 5 77 <. (all the Greek indices count beads, 
not spatial dimensions). In general, it docs not coincide 
with the conventional scalar product. If now some force 
is applied to certain bead in the packing, the above pro- 
jective procedure allows one to determine the response of 
the system, i.e. the change of the contact forces between 
all the beads below the given one. In other words one can 
follow how the perturbation propagates downward. Since 
the equations of mechanical equilibrium are linear, and 
beads are assumed to be rigid enough to preserve their 
sizes, the response of the system to the applied force is 



also linear. This linearity can be destroyed only by vi- 
olating the condition of positivity of the contact forces, 
which implies the rearrangement of the packing. While 
the topology (and geometry) of the network is preserved, 
one can introduce the Green function to describe the re- 
sponse of the system to the applied forces. Namely, force 
fx applied to certain bead A results in the following addi- 
tional force acting on another bead, fi (lying below A): 

f M = G mA • f A (12) 

Here Gx^ is a tensor Green function, which can be 
calculated as the superposition of all the projection se- 
quences (i.e. trajectories), which lead from A to fi. 

The stress field fj y in the system of frictionless spher- 
ical beads can be introduced in the following way [ p7| |: 

o- y 00 = ^ H /«/3</3 n i/3 i? «/3 (5 ( x « ~ x ) ( 13 ) 

Here R a p — \x a — xp\. As we have just shown, the 
magnitude of the force f a @ transmitted along the con- 
tact unit vector n Q( g can be expressed as an appropriate 
projection of the total force F a acting on the bead a 
from above. This allows one to express the stress tensor 
in terms of the vector field F a : 

a%J ( x ) = 53 ( F ^\ n a f3) a n l a/3 n : > a(3 R a f 3 S(x a - x) 

(14) 

In order to obtain the continuous macroscopic descrip- 
tion of the system, one has to perform the averaging of 
the stress field over a region much larger than a bead. 
At this stage we make a mean-field approximation for 
the force F Q acting on a given bead from above: we re- 
place F a by its average F over the region. To be valid, 
this assumption requires that 

a/3 

« £ (*W} Q <P<p R «P (15) 

a/3 

For certain simple geometries, the mean-field approxi- 
mation is exact. One example is the simple square lattice 
treated in Section [lj[ In any regular lattice with one bead 
per unit cell, all the F a 's must be equal under any uni- 
form applied stress. Thus replacing F Q by its average 
changes nothing. If this lattice is distorted by displacing 
its soft modes, the F a are no longer equal and the validity 
of the mean-field approximation can be tested. Figure ^ 
shows a periodic distortion with fourbeads per unit cell. 
For example, under an applied vertical force, the bottom 
forces oscillate to the left and right. Nevertheless, the 
stress crossing the bottom row, like that crossing the row 
above it, is the average force times the length. One may 
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verify that the F Q may also be replaced by its average 
when the applied force is horizontal. Though the mean- 
field approximation is exact in these cases, it is clearly 
not exact in all. In the lattice of Figure || the mean field 
approximation may be inexact if one considers a region 
not equal to a whole number of unit cells. 

A disordered packing may be viewed as a superposition 
of periodic soft modes like those of Figure ||. Each such 
mode produces fluctuating forces, like those of the exam- 
ple. But after averaging over an integer number of unit 
cells, the stress may depend on only the average force F. 
A disordered packing need not have a fixed co-ordination 
number as our example does. This is another possible 
source of departure from the mean- field result. 

Now, it becomes an easy task to perform a local aver- 
aging of the Eq. |l4| for the stress field in the vicinity of 
a given point x, replacing F Q by its average: 



cr'y'(x) = pF fe (x)r fcy '(x) 



(16) 



Here p is the bead density, F(x) is the force F a averaged 
over the beads a in the vicinity of the point x, and the 
third-order tensor f characterizes the local geometry of 
the packing: 



Mi 



(x) 



la/3; 



(17) 



This equation is similar in spirit to one derived by Ed- 
wards for the case of a d + 1 co-ordinated packing of 
spheres with friction p3|. Our geometric tensor r plays 
a role analogous to that of the fabric tensor in that treat- 
ment. 

The stress field satisfies the force balance equation, Eq. 
(||). Since this is a vector equation, it normally fails to 
give a complete description of the tensor stress field. In 
our case, however, the stress field has been expressed in 
terms of the vector field F. This creates a necessary clo- 
sure for the force balance equation. It is important to 
note that the proposed macroscopic formalism is com- 
plete for a system of arbitrary dimensionality: there is a 
single vector equation and a single vector variable. We 
now discuss the application of the above macroscopic for- 
malism in two special cases. First we consider the equa- 
tions of stress propagation in two dimensions. Then we 
discuss a packing of arbitrary dimensionality but with 
uniaxial symmetry. It is assumed to have no preferred 
direction other than that of gravity. 



A. Two-dimensional packing. 

In two dimensions, according to Eq. (|l6|), the stress 
tensor a can be written as a linear combination of two r 
tensors. 



(7 = Fl<7l + F 2 a 2 , 



(18) 



where [cri]^ = r ly and [o^] 13 = t 2%3 . Since the <j\ and a 2 
are properties of the medium and are presumed known, 



the problem of finding the stress profile a(x) becomes 
that of finding F\ and F 2 under a given external load. 
Rather than determining these F's directly, we may view 
Eq. (18) as a constraint on a. The form ( |18| ) constrains 
a to lie in a subspace of the three-dimensional space of 
stress components a = (a xx , cr yy , a xy ). It must lie in the 
two-dimensional subspace spanned by a\ and a 2 . This 
constraint amounts to one linear constraint on the com- 
ponents of cr, of the form 







(19) 



where the u tensor is determined by &i, and a 2 . Specifi- 
cally, u may be found by observing that the determinant 
of the vectors <r, o\ , a 2 must vanish. Expanding the de- 
terminant by minors to obtain the coefficients of the cr lJ , 
one finds 
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(20) 



Eq. |l9J has the same "null-stress" form as that intro- 
duced by Wittmer et al [||, whose original arguments 
were based on a qualitative analysis of the problem. By 
an appropriate choice of the local co-ordinates (£, n), the 
u tensor can be transformed into co-ordinates such that 
u« = u m = 0. Then the null stress condition becomes 
fjC'f = tjnt, = 0. This implies that, according to force 
balance equation (||), the non-zero diagonal components 
of the stress tensor "propagate" independently along the 
corresponding characteristics, £ = const and n = const: 



(21) 



Our microscopic approach gives an alternative founda- 
tion for the null-stress condition, Eq. (|l9|), and allows 
one to relate the tensor u in this equation to the local 
geometry of the packing. Our general formalism is not 
limited to the two-dimensional case, and in this sense, is 
a generalization of the null-stress approach. 



B. Axially-symmetric packing. 

Generally, there are two preferred directions in the se- 
quential packing: that of the gravitational force, g , and 
that of the growth surface n. In the case when these two 
directions coincide, the form of the third-order tensor f , 
Eq. ([It]), should be consistent with the axial symmetry 
associated with the single preferred direction, n. Since 
r fey is symmetric with respect to i <-» j permutation, it 
can be only a linear combination of three tensors: n fc n z n J , 
n k 8 t: > and 8 kl n^ + 5 k; >n l , for general spatial dimension d. 

Let cr y be the stress tensor in the d-dimensional space 
= 0...d— 1, and index '0' corresponds to the vertical 
direction). From the point of view of rotation around 



G 



the vertical axis the stress splits into scalar a 00 , d — 1 
dimensional vector <7 0a (a = X..A - 1) d — 1 dimensional 
tensor a ab . According to our constitutive equation [if], 
the stress should be linear in vector F, which itself splits 
into a scalar F° and a vector F a with respect to horizon- 
tal rotations. Since the material tensor r is by hypoth- 
esis axially symmetric, the only way that the "scalar" 
a 00 may depend on F is to be proportional to "scalar" 
F°. Likewise, the only way "tensor" a ab can be linear 
in F is to be proportional to 5 ab F°. Therefore, in the 
axially-symmetric case 



a ab = \6 ab <r 00 , 



(22) 

where the constant A is eg. r 011 / T 000 . This constitutive 
equation allows one to convert the force balance equation 
(||) to the following form: 



d°CT 00 + d a a a0 = f ext - d°a a0 + \d a a 00 = f° xt 



(23) 



In the case of no external force, we may take d° of the first 
equation and combine with the second to yield a wave 
equation for cr 00 . Evidently a ab , being a fixed multiple 
of cr 00 , obeys the same equation. Similar manipulation 
yields the same wave equation for a 0a and a a0 . Thus ev- 
ery component of stress satisfies the wave equation with 
vertical direction playing the role of time and VA being 
the propagation velocity. 



V. DISCUSSION 

In this section we consider how well our model should 
describe real systems of rigid, packed units. As stated 
above, our model is most relevant for emulsions or dense 
colloidal suspensions, whose elementary units ae well de- 
scribed as frictionless spheres. Under very weak com- 
pression the forces between such units match our model 
assumptions. However, our artificial procedure of se- 
quential packing bears no obvious resemblance to the 
arrangements in real suspensions. We argue below that 
our model may well have value even when the packing is 
not sequential. More broadly we may consider the con- 
nection between our frictionless model and real granular 
materials with friction. The qualitative effect of adding 
friction to our sequential packing is to add constraints so 
that the network of contacts is no longer minimally con- 
nected. Thus the basis for a null-stress description of the 
force transmission is compromised. We argue below that 
friction should cause forces to propagate as in an elastic 
medium, not via null-stress behavior. 



A. Sequential packing 

We first consider the consequences of our sequential 
packing assumption. One consequence is that each bead 



has exactly d supporting contacts. These lead succes- 
sively to earlier particles, forming a treelike connectiv- 
ity from supported beads to supporters. Although the 
counting arguments of Section [ij show that the propa- 
gating stress approach should be applicable to a wide 
class of frictionless systems, the continuum description 
of Section IV depends strongly on the assumed sequen- 
tial order. Now, most packings are not sequential, and 
even when beads are deposited in sequence, they may un- 
dergo rearrangements that alter the network of connec- 
tions. However, it is possible to modify our arguments 
to take account of such re-arrangements. Our reason- 
ing depends on the existence of d supporting contacts for 
each bead. Further, every sequence of supporting con- 
tacts starting at a given bead must reach the boundary 
of the system without passing through the starting bead: 
there must be no closed loops in the sequence. 

Even in a non-sequential packing we may define a net- 
work of supporting contacts. First we define a downward 
direction. Then, for any given bead in the pack, we de- 
fine the supporting contacts to be the d lowest contacts. 
With probability 1, each bead has at least d contacts. 
Otherwise it is not stable. Typically a supporting bead 
lies lower than the given bead. Thus the typical sequence 
of supporting contacts leads predominantly downward, 
away from the given bead, and returns only rarely to the 
original height. A return to the original bead must be 
even more rare. One may estimate the probability that 
there is a loop path of supporting contacts under sim- 
ple assumptions about the packing. As an example we 
suppose the contacts on a given bead to be randomly 
distributed amongst the 12 sites of a randomly-oriented 
close-packed lattice. We further imagine that these sites 
are chosen independently for each bead, with at least 
one below the horizontal. Then the paths are biased 
random walks with a mean steplength of .51 diameters 
and a root-mean-square steplength of about 1.2 times the 
mean. The probability of a net upward displacement of 
1 or more diameter is about one percent. It appears that 
our neglect of loop paths is not unreasonable. 



B. Friction 

The introduction of friction strongly affects most of our 
arguments. Friction creates transverse as well as normal 
forces at the contacts. The problem is to determine posi- 
tions and orientations of the beads that lead to balanced 
forces and torques on each. If the contact network is min- 
imally connected, the forces can be determined without 
reference to deformations of the particle. But if the net- 
work has additional constraints, it is impossible to sat- 
isfy these without considering deformation. This is no 
less true if the beads are presumed very rigid. We first 
give an example to show that in a generic packing the 
deformability alters the force distribution substantially. 
We then give a prescription for defining the deformation 
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and hence the contact forces unambiguously. 

In our example we imagine a two-dimensional sequen- 
tial packing and focus on a particular bead, labeled 0, 
as pictured in Figure ||. We presume that the beads 
are deposited gently, so that each contact forms with- 
out tangential force. Thus when the bead is deposited, 
it is minimally connected: its weight determines the two 
(normal) supporting forces, labeled 1 and 2. Thenceforth 
no slipping is allowed at the contact. Later during the 
construction of the pack bead must support the force 
from some subsequent bead. This new force is normal, 
since it too arises from an initial contact. But the new 
force creates tangential forces on the supporting contacts 
1 and 2. To gauge their magnitude, we first suppose that 
there is no friction at contacts 1 and 2, while the sup- 
porting beads remain immobile. Then the added force 
F leads to a compression. We denote the compression of 
the contact 1 as S. With no friction, the contact 2 would 
undergo a slipping displacement by an amount of order 
6. Friction forbids this slipping and decrees deformation 
of the contact instead. The original displacement there 
would create an elastic restoring force of the same order 
as the original F. Thus the imposition of friction creates 
new forces whose strength is comparable to those with- 
out friction. The frictional forces are not negligible, even 
if the beads are rigid. Increasing the rigidity lessens the 
displacements 5 associated with the given force F, but 
it does not alter the ratio of frictional to normal forces. 
Neither are the frictional forces large compared to the 
normal forces. Thus a coefficient of friction /i of order 
unity should be sufficient to generate enough frictional 
force to prevent slipping of a substantial fraction of the 
contacts. 

The contact forces T\ and T2 cannot be determined by 
force balance alone, as they could in the frictionless case. 
Now the actual contact forces are those which minimize 
the elastic energy of deformation near the two contacts. 
This argument holds not just for spheres but for general 
rounded objects. 

Though the new tangential forces complicate the de- 
termination of the forces, the determination need not be 
ambiguous. We illustrate this point for a sequential pack- 
ing on a bumpy surface with perfect friction. We choose 
a placement of the successive beads so that no contact 
re-arrangements occur. If only a few beads have been 
deposited in the container, the forces are clearly well 
determined. Further, if the forces are presumed well- 
determined up to the Mth bead, they remain so upon 
addition of the M + 1st bead. We presume as before 
that the new bead exerts only normal forces on its im- 
mediate supporters. Each supporter thus experiences a 
well-defined force, as shown in Section ||. But by hypoth- 
esis, these supporting beads are part of a well-connected, 
solid object, whose contacts may be regarded as fastened 
together. Thus the displacem ents and rotations of each 
bead are a well-defined function of any small applied load. 
Once the M + l'st bead has been added, its supporting 
contacts also support tangential force, so that it responds 



to future loads as part of the elastic body. 

We conclude that a sequential packing with perfect 
friction, under conditions that prevent contact rearrange- 
ments, transmits forces like a solid body. Small changes 
in stress 5a lJ in a region give rise to proportional changes 
in the strain £7 . This proportionality is summarized by 
an elasticity tensor K ljke : da 13 — K tjke Sj ke . The elastic 
tensor K should depend in general on how the pack was 
formed; thus it may well be anisotropic. 

This elastic picture is compromised when the limita- 
tions of friction are taken into account. As new beads 
are added, underlying contacts such as contacts 1 and 
2 of Figure [| may slip if the tangential force becomes 
too great. Each slipping contact relaxes so as to satisfy 
a fixed relationship between its normal force N and its 
tangential force T: viz.\T\ = n\N\. If /z were very small, 
virtually all the contacts would slip until their tangen- 
tial force were nearly zero. Then the amount of stress 
associated with the redundant constraints must become 
small and the corresponding elastic moduli must become 
weak. Moreover, as /j. approaches 0, the material on any 
given scale must become difficult to distinguish from a 
frictionless material with unidirectional stress propaga- 
tion. Still, redundant constraints remain on the average 
and thus the ultimate behavior at large length scales (for 
a given fj,) must be elastic, provided the material remains 
homogeneous. 

C. Force-generated contacts 

Throughout the discussion of frictionless packs we have 
ignored geometric configurations with probability zero, 
such as beads with redundant contacts. Such contacts 
occur in a close-packed lattice of identical disks, for ex- 
ample. Though such configurations are arbitrarily rare 
in principle, they may nevertheless play a role in real 
bead packs. Real bead packs have a finite compress- 
ibility; compression of beads can create redundant con- 
tacts. Thus for example a close-packed lattice of iden- 
tical spheres has six contacts per bead, but if there is 
a slight variability in size, the number of contacts drops 
to four. The remaining two beads adjacent to a given 
bead do not quite touch. These remaining beads can 
be made to touch again if sufficient compressive stress 
is applied. Such stress-induced redundant contacts must 
occur in a real bead with some nonzero density under 
any nonzero load. These extra contacts serve to stabilize 
the pack, removing the indeterminate forces discussed in 
Section D. To estimate the importance of this effect, we 
consider a large bead pack compressed by a small fac- 
tor 7. This overall strain compresses a typical contact 
by a factor of order 7 as well. The number of new con- 
tacts may be inferred from the pair correlation function 
g(r). Data on this g(r) is available for some computer- 
generated sequential packings of identical spheres of ra- 
dius R ply . These data show that g(r) has a finite value 
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near 1 at r = 2R. Thus the number of additional con- 
tacts per bead that form under a slight compression by 
an amount Sr is given by SZ = 6<ftg(2R)6r / R ~ A-y. Here 
4> ~ .6 is the volume fraction of the beads. These ex- 
tra contacts impose constraints that reduce the number 
of undetermined boundary forces in a compact region 
containing M' beads and roughly M' 2 / 3 surface beads. 
The remaining number of undetermined boundary forces 
now averages l-ZVext — M'SZ. The first term is of order 
A/' 2 / 3 , and must thus be smaller than the second term 
for Af' 1 / 3 ~ {SZ)~ X . For M' larger than this amount, 
there are no further undetermined forces and the region 
becomes mechanically stable. Moukarzel reaches a 
similar conclusion by a somewhat different route. 

If the pack is compressed by a factor of 7, stability oc- 
curs for M' 1 / 3 > I/7 — a region roughly I/7 beads across. 
In a typical experiment [|29| the contact compression 7 is 
1CP 4 or less, and the system size is far smaller than 10 4 
beads. Thus compression- induced stability should be a 
minor effect here. Still, this compression-induced stabil- 
ity might well play a significant role for large and highly 
compressed bead packs such as compressed emulsions 
[^1 . In some of the large packs of Ref. ||, compression- 
induced stability may also be important. 



D. Experimental evidence 

We have argued above that undeformed, frictionless 
beads should show unidirectional, propagating forces 
while beads with friction should show elastic spreading 
of forces. The most direct test of these contrasting be- 
haviors is to measure the response to a local perturbing 
force H] . Thus eg. if the pile of Figure [y is a null- 
stress medium, the local perturbing force should propa- 
gate along a light cone and should thus be concentrated 
in a ring- like region at the bottom |jL5| . By contrast, if 
the pile is an elastic medium the perturbing force should 
be spread in a broad pattern at the bottom, with a max- 
imum directly under the applied force. Existing experi- 
mental information seems inadequate to test either pre- 
diction, but experiments to measure such responses are 
in progress PQ] . 

As noted above, emulsions and colloids are good re- 
alizations of the frictionless case. The contacts in such 
systems are typically organized by hydrostatic pressure 
or by flow, rather than by gravity. Still, our general argu- 
ments for unidirectional propagation should apply. Ex- 
tensive mechanical measurements of these systems have 
been made [p6L^4| . The shear modulus study of Weitz 
and Mason J26[ illustrates the issues. The study spans 
the range from liquid-like behavior at low volume frac- 
tions to solid-like behavior at high volume fractions. In 
between these two regimes should lie a state where the 
emulsion droplets are well connected but little deformed. 
The emulsion in this state should show unidirectional 
force transmission. It is not clear how this should affect 



the measured apparent moduli. 

Other indirect information about force propagation 
comes from the load distribution of a granular pack on 
its container, such as the celebrated central dip under a 
conical heap of sand ||. These data amply show that 
the mechanical properties of a pack depend on how it 
was constructed. Theories postulating null-stress behav- 
ior have successfully explained these data ||. But con- 
ventional elasto-plastic theories have also proved capable 
of producing a central dip Q . An anisotropic elastic ten- 
sor may also be capable of explaining the central dip. 

Another source of information is the statistical distri- 
bution of individual contact forces within the pack or 
at its boundaries. The measured forces become expo- 
nentially rare for strong forces fl7| , ^o|| . Such exponen- 
tial falloff is predicted by Coppersmith's "q model" [jl7| , 
which postulates unidirectional force propagation. Still, 
it is not clear whether this exponential falloff is a distin- 
guishing feature of unidirectional propagation. A disor- 
dered elastic material might well show a similar exponen- 
tial distribution. 

Computer simulations should also be able to test our 
predictions. Recent simulations |31 ^ have focussed on 
stress-induced restructuring of the force-bearing contact 
network. We are not aware of a simulation study of the 
transmission of a local perturbing force. Such a pertur- 
bation study seems quite feasible and would be a valu- 
able test. We have performed a simple simulation to test 
the mean-field description of stress in frictionless packs. 
Preliminary results agree well with the predictions. An 
account of our simulations will be published separately. 



VI. CONCLUSION 

In this study we have aimed to understand how force 
is transmitted in granular media, whether via clastic re- 
sponse or via unidirectional propagation. We have identi- 
fied a class of disordered systems that ought to show uni- 
directional propagation. Namely, we have shown that in 
a general case a system of frictionless rigid particles must 
be isostatic, or minimally connected. That is, all the 
inter-particle forces can in principle be determined from 
the force balance equations. This contrasts with stati- 
cally undetermined, elastic systems, in which the forces 
cannot be determined without self-consistently finding 
the displacements induced by those forces. Our general 
equation-counting arguments suggest that isostatic prop- 
erty of the frictionless packing results in the unidirec- 
tional propagation of the macroscopic stress. 

We were able to demonstrate this unidirectional prop- 
agation explicitly by specializing to the case of sequential 
packing. Here the stress obeys a condition of the previ- 
ously postulated null-stress form Q ; our system provides 
a microscopic justification for the null-stress postulate. 
Further, we could determine the numerical coefficients 
entering the null-stress law from statistical averages of 







the contact angles by using a mean field hypothesis (de- 
coupling Anzatz). We have devised a numerical simu- 
lation to test the adequacy of the sequential packing as- 
sumption and the mean-field hypothesis. The results will 
be reported elsewhere. 

If we add friction in order to describe macroscopic 
granular packs more accurately, the packing of rigid par- 
ticles no longer needs to be isostatic, and the system is 
expected to revert to elastic behavior. This elasticity 
does not arise from softness of the beads or from a pecu- 
liar choice of contact network. It arises because contacts 
that provide only minimal constraints when created can 
provide redundant constraints upon further loading. 

We expect our formalism to be useful in understanding 
experimental granular systems. It is most applicable to 
dense colloidal suspensions, where static friction is nor- 
mally negligible. Here we expect null-stress behavior to 
emerge at scales large enough that the suspension may be 
considered uniform. We further expect that our mean- 
field methods will be useful in estimating the coefficients 
in the null-stress laws. In macroscopic granular packs 
our formalism is less applicable because these packs have 
friction. Still, this friction may be small enough in many 
situations that our picture remains useful. Then our mi- 
croscopic justification may account for the practical suc- 
cess of the null-stress postulate M for these systems. 
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FIGURES 




FIG. 1. Contrast between elastic (elliptical) and unidirec- 
tional (hyperbolic) force propagation in a sand-pile shaped 
object. A: elastic response to an imposed local downward 
force. Pictured lines of force represent the current of vertical 
force, i.e., the stress contracted into a downward unit vector. 
Near the source, this field is symmetric about a horizontal 
plane through the object J^|; part of the force is transmitted 
through points above the source. B: unidirectional response 
to an imposed local downward force. The imposed force is 
transmitted to neighbors below the source, and is further 
transmitted to neighbors below these. No force is transmitted 
to points above the source. 
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FIG. 2. A sequential packing in a channel. When left wall 
is displaced outward by a small amount 5, the beads shift to 
the positions shown by dashed lines, and the top of the pile 
shifts by an amount e. 




FIG. 3. A buckled square lattice illustrating the prop- 
agation of inhomogeneous forces. Bottom row of sites has 
alternating wide and narrow spacing. Arrows indicate the 
unequal forces on these sites. 



FIG. 4. The effect of friction on a triad of beads. In the 
absence of friction, the applied force F is transmitted entirely 
to contact 1, causing a displacement S. This would result in 
a sliding displacement of contact 2 by an amount 8. With 
friction, contact 2 cannot slide; it must deform the contact 
region by an amount of order 8. Thus the applied force F is 
shared between contacts 1 and 2. The force is distributed so 
as to minimize the total elastic energy at contacts 1 and 2. 
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